home *** CD-ROM | disk | FTP | other *** search
/ Celestin Apprentice 7 / Apprentice-Release7.iso / Source Code / Libraries / VideoToolbox 97.08.16 / VideoToolboxSources / DrawGrating.c < prev    next >
Encoding:
Text File  |  1996-08-25  |  5.5 KB  |  201 lines  |  [TEXT/CWIE]
  1. /*
  2. DrawGrating.c
  3.  
  4. Quickly draws a grating (cos, sin, square, or missing fundamental) with a
  5. gaussian envelope in the current port, bounded by the supplied rect. The routine
  6. is accurate and takes advantage of even and odd symmetries for speed. It is not
  7. entirely general, e.g. orientation is restricted to vertical and horizontal.
  8.  
  9. Squarewaves have an undefined value at zero crossings, e.g. at the origin, so I've
  10. shifted all the functions by a half pixel up and to the left.
  11.  
  12. HISTORY:
  13. 11/2/94 dgp shrink by s->noiseCheckSize, since it will later be expanded by that amount.
  14. 3/26/95 dgp as requested by Manoj, there are now separate space constants along and across
  15.             bars.
  16. 6/16/95 dgp replace "noiseCheckSize" by "signalCheckSize".
  17. 6/17/95 dgp changed "signalBounds" in first line to "stimulusBounds"
  18. 6/17/95 dgp eliminated all reference to s->signalContrast, since that applies to the 
  19.             entire sequence of images.
  20. 8/25/96    dgp    eliminated need for Numerical Recipes.
  21. */
  22. //#include "TextInNoise.h"
  23. static double *vector(long low, long high);
  24. static unsigned long *lvector(long low, long high);
  25. static void free_vector(double *v, long low, long high);
  26. static void free_lvector(unsigned long *v, long low, long high);
  27. enum{kCosinewave=0,kSinewave,kSquarewave,kMissingFundamental};// gratingType
  28. typedef struct{
  29.     double degreesPerPixel;
  30.     double spatialFrequency,spaceConstantAlongBars,spaceConstantAcrossBars;
  31.     Rect bounds;
  32.     Point origin;    // origin of the functions will be half a pixel up and to the left of this.
  33.     int gratingType;
  34.     int checkSize;
  35.     Boolean vertical;
  36.     int eAmplitude,eOffset;
  37. }GratingSpec;
  38. void DrawGrating(GratingSpec *g);
  39.  
  40. #undef MAX
  41. #undef MIN
  42. #define MAX(a,b) ((a) > (b) ? (a) : (b))
  43. #define MIN(a,b) ((a) < (b) ? (a) : (b))
  44.  
  45. void DrawGrating(GratingSpec *g)
  46. {
  47.     unsigned long *gratingPix;
  48.     double *fX,*fXEnvelope,*fY,*fTemp,w,e1,b,dc,a;
  49.     int iMax,width,fXSymmetry=0,fYSymmetry=0;
  50.     enum{kEven=1,kOdd};
  51.     Rect r;
  52.     int i,j;
  53.     int leftEdge,verticalOrigin;
  54.  
  55.     r=g->bounds;
  56.     OffsetRect(&r,-g->origin.h,-g->origin.v);
  57.     ShrinkRect(&r,g->checkSize,g->checkSize);    // Shrink
  58.     iMax=MAX(r.right,-r.left)-1;
  59.     iMax=MAX(iMax,MAX(r.bottom,-r.top)-1);
  60.     fX=vector(-iMax-1,iMax);
  61.     fXEnvelope=vector(-iMax-1,iMax);
  62.     fY=vector(-iMax-1,iMax);
  63.     width=r.right-r.left;
  64.     b=1/(g->spaceConstantAlongBars/g->degreesPerPixel);
  65.     b*=g->checkSize;    // Shrink
  66.     fYSymmetry=kEven;
  67.     for(i=0;i<=iMax;i++){
  68.         a=(i+0.5)*b;
  69.         fY[-i-1]=fY[i]=exp(-a*a);
  70.     }
  71.     if(g->spaceConstantAlongBars!=g->spaceConstantAcrossBars){
  72.         b=1/(g->spaceConstantAcrossBars/g->degreesPerPixel);
  73.         b*=g->checkSize;    // Shrink
  74.         for(i=0;i<=iMax;i++){
  75.             a=(i+0.5)*b;
  76.             fXEnvelope[-i-1]=fXEnvelope[i]=exp(-a*a);
  77.         }
  78.     }else{
  79.         for(i=0;i<=iMax;i++) fXEnvelope[-i-1]=fXEnvelope[i]=fY[i];
  80.     }        
  81.     assert(fYSymmetry==kEven);
  82.     w=2*PI*g->spatialFrequency*g->degreesPerPixel;
  83.     w*=g->checkSize;    // Shrink
  84.     e1=g->eAmplitude;
  85.     switch(g->gratingType){
  86.     case kCosinewave:
  87.         fXSymmetry=kEven;
  88.         for(i=0;i<=iMax;i++){
  89.             fX[-i-1]=fX[i]=e1*fXEnvelope[i]*cos((i+0.5)*w);
  90.         }
  91.         break;
  92.     case kSinewave:
  93.         fXSymmetry=kOdd;
  94.         for(i=0;i<=iMax;i++){
  95.             fX[i]=e1*fXEnvelope[i]*sin((i+0.5)*w);
  96.             fX[-i-1]=-fX[i];
  97.         }
  98.         break;
  99.     case kSquarewave:
  100.         fXSymmetry=kOdd;
  101.         for(i=0;i<=iMax;i++){
  102.             if(sin((i+0.5)*w)>=0)fX[i]=1;
  103.             else fX[i]=-1;
  104.             fX[i]*=e1*fXEnvelope[i];
  105.             fX[-i-1]=-fX[i];
  106.         }
  107.         break;
  108.     case kMissingFundamental:
  109.         fXSymmetry=kOdd;
  110.         for(i=0;i<=iMax;i++){
  111.             a=sin((i+0.5)*w);
  112.             fX[i]=(a>0)?1:-1;
  113.             fX[i]-=(4/PI)*a;
  114.             fX[i]*=e1*fXEnvelope[i];
  115.             fX[-i-1]=-fX[i];
  116.         }
  117.         break;
  118.     }
  119.     if(!g->vertical){
  120.         fTemp=fX;
  121.         i=fXSymmetry;
  122.         fX=fY;
  123.         fXSymmetry=fYSymmetry;
  124.         fY=fTemp;
  125.         fYSymmetry=i;
  126.     }
  127.     gratingPix=lvector(-iMax-1,iMax);
  128.     dc=g->eOffset+0.5;
  129.     leftEdge=r.left+g->origin.h/g->checkSize;
  130.     verticalOrigin=g->origin.v/g->checkSize;
  131.     for(j=iMax;j>=0;j--){
  132.         if(j>=r.bottom && -j-1<r.top)continue;
  133.         switch(fXSymmetry){
  134.         case kEven:
  135.             for(i=iMax;i>=0;i--)gratingPix[-i-1]=gratingPix[i]=dc+fY[j]*fX[i];
  136.             break;
  137.         case kOdd:
  138.             for(i=iMax;i>=0;i--){
  139.                 a=fY[j]*fX[i];
  140.                 gratingPix[i]=dc+a;
  141.                 gratingPix[-i-1]=dc-a;
  142.             }
  143.             break;
  144.         default:
  145.             assert(0);
  146.         }
  147.         if(j<r.bottom)SetPixelsQuickly(leftEdge,j+verticalOrigin,&gratingPix[r.left],width);
  148.         if(-j-1>=r.top){
  149.             switch(fYSymmetry){
  150.             case kEven:
  151.                 break;
  152.             case kOdd:
  153.                 for(i=iMax;i>=-iMax-1;i--)gratingPix[i]=dc-(gratingPix[i]-dc);
  154.                 break;
  155.             default:
  156.                 assert(0);
  157.             }
  158.             SetPixelsQuickly(leftEdge,-j-1+verticalOrigin,&gratingPix[r.left],width);
  159.         }
  160.     }
  161.     free_lvector(gratingPix,-iMax-1,iMax);
  162.     free_vector(fX,-iMax-1,iMax);
  163.     free_vector(fXEnvelope,-iMax-1,iMax);
  164.     free_vector(fY,-iMax-1,iMax);
  165. }
  166.  
  167. // the following are based on similar routines in nrutil.c of Numerical Recipes in C.
  168.  
  169. static double *vector(long low, long high)
  170. /* allocate a double vector with subscript range v[low..high] */
  171. {
  172.     double *v;
  173.  
  174.     v=(double *)malloc((high-low+1+1)*sizeof(double));
  175.     if (v==NULL) PrintfExit("allocation failure in vector()");
  176.     return v-low+1;
  177. }
  178. static unsigned long *lvector(long low, long high)
  179. /* allocate an unsigned long vector with subscript range v[low..high] */
  180. {
  181.     unsigned long *v;
  182.  
  183.     v=(unsigned long *)malloc((high-low+1+1)*sizeof(long));
  184.     if (v==NULL) PrintfExit("allocation failure in lvector()");
  185.     return v-low+1;
  186. }
  187. static void free_vector(double *v, long low, long high)
  188. /* free a double vector allocated with vector() */
  189. {
  190.     high;    // "use" argument to make compiler happy
  191.     free(v+low-1);
  192. }
  193. static void free_lvector(unsigned long *v, long low, long high)
  194. /* free an unsigned long vector allocated with lvector() */
  195. {
  196.     high;    // "use" argument to make compiler happy
  197.     free(v+low-1);
  198. }
  199.  
  200.  
  201.